Coexistence of superfiuid and Mott phases of lattice bosons 
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Recent experiments on strongly-interacting bosons in optical lattices have revealed the co- 
existence of spatially-separated Mott-insulating and number-fluctuating phases. The description of 
this inhomogeneous situation is the topic of this Letter. We establish that the number-fluctuating 
phase forms a superfiuid trapped between the Mott-insulating regions and derive the associated 
collective mode structure. We discuss the interlayer's crossover between two- and three-dimensional 
behavior as a function of the lattice parameters and estimate the critical temperatures for the 
transition of the superfiuid phase to a normal phase. 
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Dilute gases of ultra-cold bosons on a lattice present 
a model system for exploring quantum phases of matter. 
Experiments in optical lattice traps have demonstrated 
controlled tunability through a quantum phase transition 
between a Mott-insulating phase which has fixed parti- 
cle number on each lattice site and a phase exhibiting 
number fluctuations [3, 0]. As has been recently ob- 
served in radially symmetric traps [H, 0| , for sufficiently 
deep optical lattice potentials the system arranges itself 
into a "wedding cake structure" in which Mott-insulating 
phases of bosons commensurate with the lattice alternate 
with intcrlayers of incommensurate bosons with fluctuat- 
ing site occupancy [!, d, 0, @| . Various questions concern- 
ing such inhomogeneous systems have so far been unan- 
swered in both theory and experiment: is the interlayer 
associated with number fluctuations a condensate? If so, 
in what temperature range is the condensate robust? Is 
there collective behavior in the interlayer, analogous to 
that seen in bulk superfluids? How does the system cross- 
over from three-dimensional to two-dimensional behavior 
as the interlayer thickness is varied? Understanding these 
issues would also be important for related avenues in cold 
atomic physics, such as the interplay of spatial inhomo- 
geneity and quantum criticality [9| , realization of robust 
states for quantum computation llOII . and the physics of 
interacting fermions on a lattice [111 ]. 

It has been established that dilute gases of bosons 
in optical lattice potentials are well-described by the 
Bose-Hubbard Hamiltonian [H, [12] in which the bosons' 
movement between sites is characterized by a tunnel- 
ing term TLj = —JYluj) a l a j related to the overlap of 
the single-particle wave functions between neighboring 
sites i and j, and the on-site interaction is modeled by 
TLu = (£7/2) n i( n i — !)■ It is the external trapping 
potential, V(ri), which is responsible for breaking the 
uniformity of the system and promoting spatial coexis- 
tence of the Mott insulating and superfiuid phases at 
large interaction [1, 0, 0, S] • Analytical treatments of the 
inhomogeneous system are complicated by the fact that 
no simple approximation of the Hamiltonian can faith- 



fully describe the entire phase space. The Bogoliubov 
approximation 13, 3] captures the condensed phase for 



large J/U but breaks down close to the Mott regions. 

The decoupled-site approximation @, m, m, eh, valid 

when J <C £7 and the boson density is close to a com- 
mensurate value, works well within and close to the Mott 
regions but fails deep within the incommensurate phase. 
A third possibility is presented by the "pseudo-spin" ap- 
proximations [It], EM 03 , valid for intermediate values of 
J/£7, which bypass these shortcomings by treating kinetic 
energy and interactions on comparable footing. 

In this Letter, we employ a pseudo-spin approximation 
of the Bose-Hubbard model to describe the inhomoge- 
neous systems where the density of bosons varies as a 
result of a confining trap. Concentrating on a single in- 
terlayer trapped between two Mott-insulating phases, we 
show that number fluctuations give rise to a condensate 
with a well-defined order parameter and derive the dy- 
namical equations governing the system. We obtain the 
collective excitation spectrum of the interlayer conden- 
sate and show that in the homogeneous limit, it properly 
reproduces the known properties of bulk superfluids. We 
explore the behavior of the collective modes as a func- 
tion of the thickness of the interlayer and show that they 
provide a signature of dimensional cross-over in the con- 
densate, which can be achieved by tuning experimental 
parameters. We conclude with a brief discussion of the 
expected mean-field critical temperature T c of an inter- 
layer superfiuid and its behavior as a function of inter- 
layer thickness. 

Focusing our attention on the Mott phases with integer 
boson filling n and n+l and the superfiuid phase at inter- 
mediate fillings, we consider a Hilbcrt space restricted to 
the number-basis states \n) and at each site. Con- 

sidering the excluded states \n — 1), and \n + 2), we find 
that their contribution to the energy is of order of J 2 /£7. 
We note that the number-fluctuations on sites are driven 
by the incommensurability of bosons with the lattice in 
the presence of the trapping potential. The truncated 
Hilbert space in the limit J/ £7 <C 1 may be represented by 
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the spin-1/2 states [l5l|l8jl, |n + l) = [ T) and \n) = \[), 
the eigenstates of the operator s z with eigenvalues ±1/2. 
The tunneling term in the Bose-Hubbard Hamiltonian 
can be identified with raising and lowering spin-1/2 op- 
erators, s + and s _ , such that a|a,- 



aj is the boson creation operator on site i. The interac- 
tion and the potential energy terms are diagonal in the 
number basis at each site and the boson number operator 
{fii = a\a,i) can be expressed in terms of the spin-1/2 ma- 
trix s z , n = n + 1/2 + s z . Thus, in the truncated Hilbcrt 
space, one obtains an effective Hamiltonian identical to 
the spin-1/2 XY model in the external "magnetic" field: 

H = -J(n + 1) £ (sf S | + s»s*) + J2(Un MiK- (i) 

(«'} ' 

Here, (ij) denotes a summation over nearest neighbors 
sites, and /ii = /i — V(y.i) defines the chemical poten- 
tial offset by the external trapping potential, V(ri). The 
chemical potential [i is set by the total number of parti- 
cles in the system, (N) = Ylii^i)- 

The pseudospin operators are coupled ferromagneti- 
cally in the x-y plane and therefore, at low temperatures 
can form an ordered state with broken U(l) symmetry 
in the plane. At the mean-field level, in the ground 
state configuration, pseudospins are aligned with the lo- 
cal "magnetic" field, B° = zJ{n + l) [2(sf ), 2(s y t ), cos 9 % ], 
where cos^ ~ (/i^ — Un)/(zJ(n + 1)), and we have as- 
sumed (si) w (sj) for nearest-neighbors. The equilibrium 
components of the pseudospin at site i are parameterized 
by angles on the sphere: 



(n + l)s+ Sj , where 



(s z ) = (1/2) cos^, (sf) = (l/2)e^sin( 



(2) 



where the angle tp independent of site index expresses 
the phase coherence in the system. The continuous de- 
generacy in the ground state is illustrated in Fig. [TJ In 
the Mott phase, the pseudospins are completely polar- 
ized along the z direction, i.e. (sf) = ±1/2, allowing the 
identification of \i± = Un ± z J(n + 1), the values of the 
chemical potential at the boundaries of the Mott states 
with n and n + 1 bosons per site (see Fig. Q}. In the xy- 
symmetry broken phase, (a^) = (s + ) / y/n + 1 ^ and we 
have a condensate with order parameter 



A= (1/Mntcr)^( 



(3) 



where Mntcr is the number of lattice sites between the 
two Mott phases. 

The locations and sizes of the interlayers can be de- 
termined by the relationship fi± = (Jt — V(r±), where 
the chemical potential, /i, is obtained self-consistently by 
fixing the total number of particles in the system, N. 
For radially-symmetric traps, a simplification occurs in 
the limit of a thin interlayer, Sr n <C r n (where r n is 
the radius at the center of the interlayer between two 
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FIG. 1: The mean-field phase diagram for the Bose-Hubbard 
Hamiltonian. The dash-dotted line corresponds to the in- 
terlayer with fluctuating site occupancy. Spontaneous sym- 
metry breaking in the ground state is shown on the sphere 
(s) 2 = 1/4: the equilibrium configuration (s) is degener- 
ate on the circle (dashed line) with nonzero order parame- 
ter |A| = (l/2)sin#. North and South poles of the sphere 
correspond to Mott states with n + 1 and n bosons per site, 
respectively. 



Mott states with particle occupation n and n + 1 and 
Sr n is its thickness). In this case, the trapping poten- 
tial can be linearized around r n and we find that the 
number of particles in the interlayer is the same as in 
the case J/U = 0, where the interlayer region would 
be filled with n and n + 1 Mott phases. Hence, the 
chemical potential at small J/U can be found by set- 
ting N — (47r/3) X)™=o 1 ( r ™/^) 3 ' wnere m is the total 
number of Mott states in the trap and £ is the lat- 
tice spacing [8j. We find that for a three-dimensional 
parabolic trapping potential, V(r) — ar 2 , the interlayer 
parameters are given by r n = (fi/a) 1 ^ 2 [1 — nU/fi] 1 ^ 2 , 
Sr n = 6J(n + l)/(ar„) when n > 0, and Sro = 3J/(axn). 
These results show that it is possible to tune the width 
of the interlayers from Sr n ~ t to Sr n ^> t, effectively 
changing the dimensionality of the layers. As a charac- 
teristic example in the range of recent experiments [H, 0] , 
a system with trap curvature a ~ h x 24Hz//im 2 , to- 
tal particle number N w 10 6 , lattice spacing 0.43 /im, 
interparticle interaction U ~ h x 10 kHz and tunneling 
strength J as h x 120 Hz hosts two Mott regions with 
n = 1 and n = 2, and two interlayers. The corresponding 
interlayer parameters are ro ~ 25 /im and Stq m 0.5 /xm; 
n ~ 14 fiva and 5r\ ~ 4 /im. We remark that for a 
harmonic trap, when gravity is taken into account, the 
system experiences a shift along the direction of the grav- 
itational field, but is otherwise unaffected. 

Having identified the interlayer region, we now turn 
to the low-energy collective modes within the inter- 
layer. These modes can be calculated using the Heisen- 
berg equations of motion for the pseudospin operators, 
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dtSi — i\H, Si]. In the mean- field approximation, one ob- 
tains the Bloch equations, 9t(sj) = (sj) xB„ where the ef- 
fective magnetic field is given by Bf = J(n+1) J2j 
(summation is over the nearest neighbors of site i), and 
Bf — zJ(n + l)cos6>i. Assuming that the character- 
istic wavelength of the excitations is much larger than 
the lattice spacing £, we approximate the sums entering 
the effective magnetic field B by their continuum limit, 



z s 



2 V 2 (s). The resulting Bloch equations 



can be analyzed as follows. The equilibrium number den- 
sity of bosons in the interlayer is po = n + 1/2 + cos 0/2 
and the number density p = po + Sp obeys the conti- 
nuity equation d t p + Vj = with the current density 
j = (J(n + l)/2) sin 2 6 Vy>. Using the relationship be- 
tween the canonically conjugate density deviation Sp and 
the phase, dt<p = —2zJ(n+ l)6p, one obtains the follow- 
ing differential equation for density fluctuations around 
equilibrium: 



<9 2 Sp = 4z (J(n + 1)£) 2 V [A 2 VSp] 



(4) 



where Ao = (1/2) sin is the local value of the or- 
der parameter vanishing at the boundaries of the Mott 
states. The form of Eq.© is identical to that governing a 
trap ped Bose-Einstein condensate in the absence of a lat- 
tice |20| with Aq playing the role of an equilibrium den- 
sity of the condensate confined between two Mott states 
to an interlayer with radius r n and width Sr n . It must 
be noted that while the equations governing density dis- 
tortions are identical to those derived from the standard 
Gross-Pitaevskii formalism [2(| for a condensate in the 
absence of a lattice, the equations of motion for the or- 
der parameter (s + ) in general do not correspond to the 
Gross-Pitaevskii form, but reproduce it in the limit of 
small density distortions. 

For the uniform case, the excitation spectrum can be 
obtained by treating the order parameter Ao as spatially- 
independent. The eigenvalue equation, Eq. ((4|, is solved 
by the Fourier transformation, Sp cx exp(ipr — iu> p t), 
where p is the wave vector. The resulting sound mode, 



cp, 



c = y/zJ(n + 1)£\ sin0|, 



(5) 



is related to the spontaneously-broken symmetry in the 
ordered state. According to the Landau criterion, the 
sound-like spectrum of Eq. §5§ makes the ordered state a 
superfluid. One notices that the speed of sound, c, goes 
to zero as one approaches the Mott phases at sin 6 = 0. 

In the trapped geometry, an estimate of the excitation 
spectrum can be obtained from the quantization condi- 
tions imposed on the wave vector p in Eq.®, with the 
speed of sound approximated by its value in the center 
of the layer, cq — yfzJ(n + 1)1. For a spherically sym- 
metric trap, the excitation modes are confined within an 
interlayer centered at radius r n with width 5r n — 2a n . 
The excitation modes terminate at the boundaries of the 
Mott regions, i.e. pj — j/a n with j = 0,1,..., which 



gives a radial mode spectrum uij ~ il r j with character- 
istic frequency Q r = J(n + l)£/a n . The quantization 
of the surface modes is related to the angular momen- 
tum L = 0, 1, . . . through p^ = L/r n which leads to 
the spectrum ujl — Q a L with characteristic frequency 
fl a = Cl r a n /r n . The degeneracy of the surface modes is 
(2L + 1) for each value of L. The perturbative calculation 
of the modes in Eq.Q in the limit a n /r n -C 1 confirms 
these estimates and gives the following result: 



j(j + 1) + 



i 



L(L + 1) 



2r? 



(2j-l)(2j 
1 



3) 



(2j-l)(2j + 3) 



(6) 



where il = \/6^r, j = 0,1,..., L = 0,1,..., and 
j' + L / 0. In the continuum approximation, the wave- 
length of the modes should be much larger than the 
lattice spacing, which sets upper bounds on the quan- 
tum numbers: L <C r n /£ and j <C a n /£. The sec- 
ond term in Eq.© is independent of L and is associ- 
ated with the curvature of the interlayer; it vanishes at 
j = 0. The lowest energy modes for thin interlayers, 
a„/r„ <C 1, correspond to angular excitations (J = 0) 
given by uj l = 2fi Q yj L{L + 1), L = 1, 2, . . .. 

We note that the mode spectrum of Eq. ((6J corre- 
sponds to that of a condensate confined with an ex- 
plicitly shell-shaped trap (for instance, a "bubble trap" 
in Ref. [2lj) since the "effective confining potential" in 
Eq. (|4]) has the form V e ff oc (r — r„) 2 /a 2 for thin in- 
terlayers. The calculation of the radial (L = 0) modes 
with j = 1,2 in Ref. [22| confirms this connection for the 
lowest-lying radial modes (analogous to "breathers" in 
spherical condensates) . 

The characteristic frequencies fl r , fl a of the radial and 
angular modes set temperature scales at which the spec- 
trum in Eq.© becomes quasiclassical, j, L 3> 1. For 
the aforementioned experimental parameters, the cor- 
responding energy scales are of the order fi r ~ 5nK, 
£l a ~ 0.5 nK. The energy of the system at finite tem- 
perature is obtained through quantization of the col- 
lective modes, E(T) = ^Z L j{2L + l)ojLjiiLj, where 
riLj = l/(exp(u>Lj/T) — 1) is the thermal occupation 
of the bosonic modes with spectrum given by Eq. ([6]), 
and the factor (2L + 1) takes into account the degener- 
acy of the angular modes. There are three temperature 
regimes in this case. In the extreme low-temperature 
limit, T <C f2 a , thermal excitations are gapped, i.e. 
E(T) cx fl a exp(—2^/2fl a /T). At intermediate temper- 
atures, Q a «T« H r , the radial modes are frozen and 
only the two-dimensional angular modes contribute to 
the energy, E(T) oc T 3 /f2 2 . At higher temperatures, 
T 3> i both radial and angular modes are excited, and 
the energy has a three-dimensional phonon-like temper- 
ature dependence, E(T) oc T 4 /f2 3 . The separation of 
temperature scales achieved by changing the interlayer 
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width from Sr n ~ £ to Sr n 3> £ tunes the effective dimen- 
sionality of the system from two to three dimensions. 

At finite temperatures, the order parameter intro- 
duced in Eq.Q is depleted from its zero temperature 
value by the collective modes. In the low-temperature 
regime, T -C J(n + 1), for wide interlayers, 5r n £, 
the order parameter depletion is similar to the case of 
a three-dimensional weakly-interacting BEC, SA(T) oc 
(T/J(n + l)) 2 . At higher temperatures, the long-range 
order is destroyed by the quasiparticle excitations whose 
wavelength is of the order of the lattice spacing. The 
critical temperature ought to be of the same order of 
magnitude as these excitations with the energies of order 
J(n + 1) (obtained from setting j ~ a n /£ and L ~ r n /£ 
in Eq.JB])). A mean-field calculation similar to the one in 
Ref. [l7| confirms the estimate and provides a mean-field 
expression for the critical temperature T^d — 3 J(n + 1). 
Thermal properties of narrow interlayers, Sr n ~ £, are 
qualitatively different. In this case the angular excita- 
tions play the dominant role. In accordance with the 
Mermin-Wagner-Hohenberg theorem [23] , the long-range 
order is destroyed but the system retains power-law cor- 
relations in the phase of the order parameter. In the 
limit that the interlayer width is comparable to the lat- 
tice spacing, Sr n ~ I, a simple model capturing the prop- 
erties of the two-dimensional system involves only phase 
variables and leads to the effective Hamiltonian, Ti v = 
(K/2) J d 2 x (V<p) 2 , where the integration is over the sur- 
face of the spherical layer, and we have introduced the 
phase stiffness K = J{n + l)/2. The Kosterlitz-Thouless 
(K-T) transition [24| between the high-temperature nor- 
mal and the low-temperature superfluid state occurs at 
temperature T 2D = (ir/2)K = (ir/4)J(n + 1). At in- 
termediate widths, Sr n > £, the phase stiffness is ap- 
proximated by K — (J(n + \)/2){5r n /£) sin 2 9 with 
sin 2 9 = (1/Sr n ) J dr sin 2 9. For the interlayer in the trap 

sin 2 9 — 2/3, and the critical temperature of the K-T 
transition is given by T c — {-k /&){8r n / £)J{n + 1), which 
is a linear function of the interlayer width, interpolating 
between two-dimensional and three-dimensional limits, 
T2D < T c < T^d- In the range of current experiments, 
for J w hx 120Hz, one obtains an estimate of the critical 
temperature, T c ~ 10 nK. 

In conclusion, we have shown that the interlayer with 
fluctuating site occupation confined between two Mott 
states becomes superfluid at low but experimentally ac- 
cessible temperatures. Employing the pseudospin model, 
we have identified the effective potential confining the su- 
perfluid and analyzed the low-energy excitations in the 
system. We have demonstrated that the effective dimen- 
sionality of the interlayer can be changed by tuning ex- 
ternal parameters. As an example of the ensuing physics 
we have suggested that the critical temperature inter- 
polates between two-dimensional and three-dimensional 



limits as one changes the width of the interlayer. A clear 
experimental signature of the interlayer condensate, ei- 
ther through time-of-flight and interference experiments, 
excitation of collective modes or radio-frequency spec- 
troscopy is yet to be obtained. 
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